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ABSTRACT 


This paper presents an accuracy analysis of a suggested approximate 
confidence interval for system maintainability Spin ee ers. 

Technically, the simulation demonstrates feasible ranges of parameter 
applicability for a fit of linear combinations of generated gamma variates 
to the gamma distribution, using the method of moments. 

The simulation has application to the classical confidence interval 
for mean time to repair of a series system, under the assumptions: of 
gamma distributed repair times, and method of moments estimators.. 

The paper provides no validated conclusions although it does display 


parameters and ranges of apparent extremely high model validity.. 
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I. INTRODUCTION 


The purpose of this paper is to present a computer simulation. in 
order to investigate the adequacy of a suggested <i confidence 
interval for system maintainability parameters. In a technical sense, 
it will demonstrate feasible ranges of parameter applicability for a fit 
of linear ccmbinations of generated gamma variates to the gamma distribu- 
tion, using the method of moments. The application of the simulation is 
to the classical confidence interval for mean time to repair a systen, 
under the assumptions of series components, gamma distributed repair 
times, and method of moments estimators. The parameters of the distribu- 
tions that were investigated were selected so that the maintainability 
curves would be representative of those curves presently of interest ir 


industry. 





IIT. SUMMARY 


A procedure for a fit of linear combinations of generated gamma 
variates to the gamma distribution, using the method of aeMemNES is 
presented in Chapter III. A simulation of this procedure is proposed 
in Chapter IV and the results are tabulated in Tables I and II of 


Chapter V. 





III. THE STATISTICAL MODEL 


A. SERIES ASSUMPTION 

Suppose we are given a system of components which are ce necessarily 
series, where the time to failure of the ith component of the system has the 
exponential Or.) distribution. 

Then, 

1) for the case where there are k, identical components in series of 
type i within the system, set re - kX ee 


2) for the case where there are 7 identical components in parallel 
| h. 
j 


Therefore, for some parallel-series combination there is obtainable 


of type j within the system, set a = i 


a series system of uniquely different comsonents whose failure ratios <ére 
assumed to be exponential. 
Thus, the system being described is a series system with k differert 


ft rt 
thm Ay See. 6M. 6UWwhere 


components and failure rates Ay do» a : j r AL 


all the \. are in the Same units. 


B. FAILURE RATE ASSUMPTION 
Let X. denote the time to failure of component i, where xX has the 


exponential Or.) distribution. Thus, 


Z Ae a 
. oh (x, ; r) de 7 2 (3-1) 


We further assume that the system fails when exactly one component 
fails or the component which caused the system failure takes the longest 


time to repair. 





C. GAMMA ASSUMPTION 
Let T; denote the time required to repair the ith component and 
suppose that T, has the gamma distribution. Thus, T, ~T (t,; a.,.b.).. 
it ai 1 ies, Sei 


Therefore, 


ee" i & 
4 “3 Ge ay T(b,) (2) 
i a 
is the density equation with 
De 
E(T,] - a, We | (a 


which is the mean time to repair (MITR) for the ith component. 


D. SYSTEM MTTR OBTAINED 
Let 9 denote the mean time to repair for the system. Therefore,,. 
6 = © P [component i fails first] ue (3-4) 
i 


and, if the assumption of a series system is valid, then 


it oo 


a 


and 
is ge 

icc) (=) (3-5) 
: r 2h 
i=l i 


© 
il 


E. CONFIDENCE LIMIT OBTAINED 


An estimated upper confidence limit for system maintainability, denoted 


a 


by 6 has been proposed. ee was derived by the following procedure: 


For any component i, coe T. is a random sample on T., 


1 ll Ah 
ele + 13” in. 
the time to repair component i. Failure data is also available on component 


aA 


i, so that an estimate d. Or r. is possible. 





aA 


Thus, 9 can be estimated by 9 where 


7 k A, 
@= ft ae = (3-6) 
i=l 
and 
io = ooo 
4 He 
Lj 


Moreover, because of the large amount of industrial testing,. it: can 


aA 


be assumed that A. = i> and for purposes of derivation we shall treat 


the r. as though they are constants. Thus, 


5 ne : 
§ = 2 ; Ts? (3~7) 
i=l 
rs Ki) by 
Efe} = 2— ~~) (3-8) 
2 AK a, 
i=l re 
and 
“ EK 2%. 2 b, 
vara ae) (—). (3-9) 
E r Ye 
i=l a. 


We shall fit a two parameter gamma to the density of 9 and obtain 
the upper eonfidence limit frem this fitted distribution. Thus we are 
assuming 98 ~ [(a,b) and use will be made of the method of moments: to make 


the fit. Through use of formulas (3-8) and (3-9) 


a k A. b. b 
Efe]= r— = (3-10) 
. aA ae. a 
i=l 1 
and 
Pa k h; a be. ty 
var (60) = = €&> ey = — (6-11) 
X 2. 2. 
i=1 ; a 


i 





The parameters a and b can be solved for, as follows: 





k r. by 
»? ren) Guw) (3-172) 
. a 
aoe i=l] i 
k a b. 
u GC oars 
i=1 a. 
i 
and 
k rs b. - 
cS Xa, 
DS ma By (3-13) 
ae 2 
i=l a, 
af 
For the available data, this then becomes: 
n Z 
KAS 
—— Te 
b = fell a 5 (3-14) 
k r. 9 
Dac) S. 
i=l A 
where, 
M4 
35 = i os (T, eee T)? 
n, j=l J 


If b is an integer, 20a has the ae distribution: therefore,. define 
(2b) as the integer closest to 2b and the following approximate confidence 


interval can be formed: 


10 





— 
i 

R 
| 


een > Bes cs)! (3-15) 


2 
20a X (1-0) (2b) 


= P Uap) (2b) | 


ae [2 < sa 


X(1-a) (2b) 


b . @ (2b) 


=P ie ] (3-16) 


Pes (2b) 
where the approximate equality compensates for the use of estimators. 
Finally, the choice of a = .20 results in the following: 
go 2p [2 < ob) _), (3-17) 
X(,80) (2b) 
Thus, the desired suggested approximate &0% confidence limit for system 


Maintainability is 8(2b)/x¢ 8) (2b) ° 


fal 





IV. THE SIMULATION PROCEDURE 


As explained in Chapter III, the system to be simulated consists of 
k components in logical series with system maintainability (8) expressed 
as: 


k i, 
oD 


b. 
eet me. 
eet GG (4-1) 


i=1 

where b./a; is the true maintainability of the ith component of the 
system. 

Denote an upper confidence limit for 9 by Q. ibe @ is in. fact. the 
exact 100 (l1-a)% upper confidence limit for 6, then 
p [2 < 8 (2b) } 2a Se 

X (1-2) (2b) 

holds to a reasonably close approximation. 


In fact, b/a should then be the ath percentile point of the simulated 


distribution of Oke 





b 
a 


as the 20th percentile point of 


The choice of a .20 defined o. 


(2) 


the distribution of 9. 


G4 





«a 


In order to investigate this distribution, a computer was utilized 
to generate the required gamma variates. Then, 500 values of 6 were 
u 


computer and ordered such that 


a “a a a 


ree eg) 500 (4-2) 


Since it was desired to display the 80% upper confidence limit for 


8 (0 > which implies 80% of the ordered values be greater than @ 


uG2 wG2) 
the 20th percentile point of the a distribution was found. This 20th 
percentile point is the 100th ordered value in formula (4-2) above and, 
if the procedure is correct, should approximately equal 2 

The primary measure of accuracy for the simulation will be the value 
of 
1.100 7 A 


b/a @-3) 


which is an expression relating the estimated value of system maintain-—- 


ability (6 ) to the gamma value for MITR (b/a). Thus, the accuracy 


ul00 
of the simulation is presented in the notion of relative error. 


Also, the statistic 


[#6 > b/al 


500 ~~ 


will be computed in order to display the relative error between the true 
value of system maintainability (b/a) and the number of generated esti- 
mates of this value (64). 

The analysis for this proposed method was conducted using different 
combinations of values of the gamma input parameters, varying . values ,. 


one of two values for k, and n, values of 10, 20, 50, and I00. 


3 





A. THE SIMULATION PROGRAM EXPLAINED 


Available in Appendix E is a flow chart for the computer simulation.. 


An explanation of the blocks on that diagram follow: 


1) Dimension the computer matrix as required to include the possible 


values for the following input parameters: 
k 


hy do» oe k 


(a, > by); (a,, bi), "7%, fa, b.) 


Nl nh eee 19 Ge 
2) iS ae 


2) Generate the following n. random repair times according to the 


gamma distribution with parameters (a. , De eetttice the (bM Scientific 


Subroutine Package GAMMA to get: 
Wine Ung See ae In 


At eT 


To? T22° 193? 2n, 


Tal? Seno See ae “kn, 
3) For each row above, compute the mean value GE) and the variance 


(3,7) EGucen tie: Dalrs: 


—~ 2 = 2 ~ 2 
ec: T,s S, S 


4) Utilize formula (3-14) to compute the method of moments estimate 


of the unknown parameter b, denoted by b. 


5) Compute the integer closest to 2b for use with the «Cob 


distribution. 


14 





a 


6) Utilize formula (3-6) to compute the estimate of system maintain- 
ability 6, denoted by Q. 

7) Utilize formula (3-16) and the ne formula (see Appendix C) to 
compute the approximate upper confidence limit on system maintainability. 
pene this value as Q 

8) Repeat step 2 through step 7 a total of 500 times. 

9) Utilize the IBM SSP SHSORT to order the 500 values of @ from 
smallest to largest. 


10) Pick out the 100th value of the above ordering (0. 0 and print 


30} 
out this value as the 20th percentile point of the distribution.. 

11) Compute the value of b/a and then utilize formula (4-3) to compute 
the primary measure of accuracy for the simulation. Print out this value. 


12) Utilize formula (4-4) to compute another measure of simulation 


accuracy. Print out this value. 


US) 





V. RESULTS 


For the procedure as presented in Chapter III, the value of system 


maintainability will follow the gamma distribution. The parameters. of 


the distribution were chosen so that the maintainability of each compon- 


ent was at a preselected value position on the curve or would be repre- 


sentative of those curves currently obtainable in industrial applications. 


(See Appendix A.) 


CASE ? 


T 


Il 


III 


IV 


n, 
al 


Nn, 
ag 


The following cases were simulated: 


INPUT PARAMETERS 


1S 
= 10, b, =1 
i 
(0) 
= 105720, 50, 
iS 
= 105-b.. = 2 
th 
= .005 
=e 105 20;,. 50, 
iS 
= 5, b, =l 
i 
= .005 
a0 20, (50), 
15 
= 30, b. = 1, 
i 
= 10, b. =1, 
i 
O05 
i 
=O... 50-5 


100 


100 


100 


100 
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CASE # INPUT PARAMETERS 


V KA 
a, = 30, bs —jiaeie= 1, -**,..10 
eee —o teed — ll, ***, k 
pe L 
i = .01 
n, = EO, 20, 50, 100 
VI k = 30 


ee oe Dy =e = 1, -**, 10 
i ft 


a, = 10, b, =1, i=11, ---, k 


i i 
A, = .005 
i 
n. = 10, 20, 50, 100 
VII k = 30 
Ae —esOembe ol st = 1. ts eye 10 
a aL 
ane — t0- eb. = te t= ll, ***, k 
af i 
A. = .Ol 
i. 
n, =7 10-720 .450% 100 
VIL k =15 
a, = 5,b.=1 
sl a 
A. = .05 
i 
n, = 10, 20, 50, 100 


Each case is specified by the relevent input parameters and the vary- 
ing number of random repair times generated (n,). The cases differ 
most significantly in their values for the parameters of the Gamma 
distribution. However, the value of rs is varied as is the value of 
k, the number of components in series. Of course, all of the cases were 


run over an identical range for n,- 
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a 


The results are compiled in the following tables. The column 
headings are self-explanatory and their use has been discussed earlier 


in this report. 
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APPENDIX A 


THE GAMMA DISTRIBUTION 
Use of the gamma function to describe maintenance-time functions. 
has been described in the literature. 
The gamma distribution is a two-parameter function. 


—-a.t, 
b. ine se 8 
a. 1 Cc. TL a3 


| r(b.) 


The parameters are denoted by the letters a and b. Of the two,.b is 
considered to be the more critical because it controls the shape of the 


curve; while a merely determines the scale of the axes. 


TO 
b=0 
.8 
.6 
b=1 
4 b=3 = 
2 
M2 3) (4 5 6 7 8 9 


Figure 1. Gamma density function; 
ao Seb = Ooo 


For the distribution as presented, the parameters a and b are 
restricted by the inequalities a> 0, b > 0. Also, both a and b will 


remain as whole numbers. 
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a 


The b parameter is seen to vary as a function of the repair concept. 
When the repair concept being described involves relatively simple 
actions (e.g., black box replacement, etc..) and very few long downtime 
tasks, the b parameters will be a small mimber. As the downtime: density 
involves more time-consuming actions, the b parameter will take oma 
larger value. 

Throughout the literature it has beem shown that the gamma: distribu- 
tion could be used to describe a family of distribution curves varying 
from the negative exponential to the normal.. Thus, by suitable choice 
of the parameters a and b, the downtime of a wide variety of reasonably 


complex equipments can be described. 
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APPENDIX B 


THE EXPONENTIAL DISTRIBUTION 
One of the most widely used distributions in fields of Reliability/ 
Maintainability is the one-parameter exponential function.. 


The function is defined by: 


£(X, , .) = r.e xX > 0 
= 0 ‘elsewhere 


where ds is the failure rate associated with component i.. 


F(X, ) 


tt 
rd 
~7-~ 
>< 
| A 
* 
me 
New! 
i 
> 


—r.t, —r\ xX. 
{ e a dt. = l-e ater Me > 0) 
0 


= 0 Gi 8, 
LS Soe 
is the cumulative distribution function, and 


E[X,] = 1/d,. 


Thus, the expected value equals the reciprocal of the parameter i: 
an 0,00 Em Wi, : 
ea i 


1.0 


£(X.) 


Figure 2. The exponential distribution 


Zo 





NOTES: 


1) The exponential distribution is a special case of both the gamma and 


weibull distributions. 


2) It is the distribution which is expected when the mechanisms are so 


complex that many deteriorations with different failure rates are 


operable. 
3) When parts have an exponential failure distribution, the equipment 


consisting of these parts also has an exponential distribution function. 
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APPENDIX C 


THE CHI-SQUARE DISTRIBUTION 
A.quite useful case of the gamma distribution occurs when a = 1/2 
and b = n/2, where n is a positive integer. Then, a one-parameter 
family of distributions is obtainable with density function: 
(n/2 - 1) -2/2 
£(2) <= = hen 
2 T (n/2) 
A random variable z having the above density is said to have the 
chi-square distribution with n degrees of freedom (denoted by er 
In Figure 3 below, the density function forn = 1, 2, and n > 2 is 


shown. 


E(z) 





Figure 3 
with E{z] =n, var (z) = 2n. 
Our interest in the chi-square distribution is based upon its many 
important applications in statistical inference. 
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The chi-square distribution is tabled for degrees of freedom up to 
30. Thus we may find in the table that value, denoted by 7 MGe and. 
a) (no 


satisfying P(Z <x.) =O, OF—8a =< 1. 


f(z) 


2 
X (a) (n) 


Figure 4. The chi-square distribution 


HoreEhe case where the degrees of freedom exceed 30 (m > 30) the 
chi-square distribution can accurately be approximated by the normal 
distribution as indicated in the following theorem: 

Teeter: Suppose that the random variable Y has distribution Xa" 
Then for sufficiently large n, the random variable /2Y has approximately 


the distribution N(V2n - 1, 1). Therefore, 


26 





BCE <tc) 


P(V2Y¥ < V2t) 


P(Y2Y - Y2n - 1 < V2t - /2n — 1) 
o(¥2t - Y2n - 1) 


where ® values are obtained from the normal tables. 
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APPENDIX D 


THE METHOD OF MOMENTS 

The oldest general method for generating estimates Bes anita 
parameters, given a sample, is known as the method of moments. It is 
a comparatively simple method and it generates "reasonable" estimators. 

In general, the method proceeds as follows: 

Given a random variable X, with distribution function FY where this 
distribution is indexed by the unknown parameter X. Assume the first 
moment of X (its mean value) is dependent in some simple way upon \, 

like Hy = g(\). Then, given a sample of n values of X, define the first 
sample moment as X. The method of moments now specifies that X be equated 
to g(\) and finally to solve for A This resulting value for . is the 
method of moments estimator of i. 

For the two dimensional case, such as X being distributed according 
to the gamma distribution with unknown parameters a and 8,. there is 
required a multi-dimensional parameter space = EB which in the gamma 
case has k = 2. Then, this space can be defined such that 


Q = {(a, 8); a > 0, B > O} G E, and further, it can be shown that. 


2 


2 2 DD. 
Hy = a8 and Uy = a8 ta®6 


where wu : is the second moment with respect to X. Then,. as in the one 


X 


dimensional case, a random sample is required which can be analyzed to 
give: 
X = af 


Sate x ee eens 


mM hO 
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which allows, as shown in Chapter III 


2 

a x2 a Sy 

a Ze) B st 
Sy x 


as acceptable estimators of the parameters a and 8. 
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APPENDIX E 





COMPUTE 
b/a 


SUBROUTINE \ | GENERATE 


AMMA GAMMA 
7 VARIATES 














COMPUTE 
8 


PRINT OUT 
SLAP ISI iCS 


ULILIZE x2 
TO COMPUTE 





ORDER 500 
VALUES OF 


SUBROUTINE 
SHSORT iy, 


Gi 


PRINT 20th 
PERCENTILE 
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1,K) 


1,K), (ALPHA(JS),BETA(JS),-J= 


E EWeEREST 


EG) 


00) ,S2(100),TU(5C0) ,ALPHA(100),BETA(190 
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DO 1234 JJ 


OOOO YOUO 


(BETA(I),ALPHA(T),X) 
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et et OE 
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GAMMAS GENERATED AND THE 
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5 
g 
g 
C 
& 
g 
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& 
6 
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12 34 
g 
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& 
g 
99 
& 
C 
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c 
C 
C 
222 
100 
195 
200 
ie 


CONTINUE 

COMPUTE VALUE FOR ESTIMATE OF PARAMETER BETA, CALLED BHAT 
BHAT=(SUM*SUM) /SUMSQ 

COMPUTE INTEGER CLOSEST TO TWO BHAT 

IH=2.0* *BHAT+O. 5 


BHAT2 = 

Serotec aot IMATE OF MEAN TIME TQ REPAIR FOR SYSTEM,.CALL IT THETA 
THETA = 0.29 . 

DO 6 I=1,K 

Meet A — THETA + (LAMBDACT OV LAMI=T S&C 1) 

CONTINUE 

See CHI S@ FORMULATION 0 Get SPATISTIC FOR CONFROENCE LIMITS 
CHISA=(2.0/(9.9 *<BHAT2)) 

CHI SQ=BHAT2*(1.C-CHISA-C. 241 7T8*SQRT(CHISA) )©*3 

TUC SJ) = CTHETA* BH AT2)/CHISQ 

CONTINUE 


NOW THAT THE 500 ESTIMATED VALUES HAVE BEEN GENERATED, SORT THEM 
INTQ ORDER 


KEY 

CALL SHSORT(TU KEY. KK} 

eee oteANOee ni Nt tHE 2ORHePERCENTILE ACCORDING TO THE FOLLOW- 
ING FORMA 

Peete (GO, 227) NiyKyLA4BDA{ LJ 

WRITE(6,200) TUC1LOO) 

meme tt it, "Ni = *15, ° st5o 9’ LA MBDA=*,F1O0~ Sy9//7/) 

meen (15/7 15F5.3 / 1555. 3 F(LOEF 4. 99F 4e0))) 

FORMAT (15) 

BU eheixg: THE 20TH PERCENTILE IS *,F15.5) 


WRITE(6,400)(TU(CT),I1=1 5,500) 
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